
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
# R code for this Project starts here
#&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
Wdistance <- read.csv("distanceW85Rev.csv")
#load("WAggregate39.RData")

load("trade011312.RData")
traderep <- trade011312
load("WmWTOlist.RData")

attach(traderep)
library(lme4)
library(foreign)
library(nlme)
#library(pcse)
library(lmtest)
library(plm)

# panel data formatting
trade <- pdata.frame(traderep, index = c("cowcode", "year"), drop.index = TRUE,row.names = TRUE)

#####################
# Type I interaction
#####################
#panel data analysis
type1.trade <- plm(newtar~l1polity2+ l1log.pop + l1log.gdppc  +l1bpc+l1ecris+l1yrsoffc+
                   l1gatt.wto, data=trade, effect="twoways", model="within")

summary(type1.trade)

# multilevel modeling counterpart
mlm.type1 <- lmer(newtar~l1polity2+ l1log.pop + l1log.gdppc  +l1bpc+l1ecris+l1yrsoffc+l1gatt.wto+(1|cowcode) + (1|year), data=traderep)
summary(mlm.type1)
                  
                  
#####################
# Type II interaction
#####################
#panel data analysis
type2.trade <- plm(newtar~l1polity2+ l1log.pop + l1log.gdppc  +l1bpc+l1ecris+l1yrsoffc+l1gatt.wto + l1usheg + l1avetar + Eight + Wmember + wtrade, data=trade,  model="within")

summary(type2.trade)

# multilevel modeling counterpart
mlm.type2 <- lmer(newtar~l1polity2+ l1log.pop + l1log.gdppc  +l1bpc+l1ecris+l1yrsoffc+l1gatt.wto + l1usheg
                   + l1avetar + Eight + Wmember + wtrade+(1|cowcode) + (1|year), data=traderep)
summary(mlm.type2)

# use AIC or BIC to compare models. or Ajusted R square

######################
# Type III Interaction
######################
#panel data analysis
type3.trade <- plm(newtar~l1polity2*Wmember+ l1log.pop + l1log.gdppc  +l1bpc+l1ecris+l1yrsoffc +l1gatt.wto + l1usheg
                   + l1avetar + Eight + wtrade, data=trade,  model="within")
summary(type3.trade)
# Multilevel models
trade.mlm <- lmer(newtar~l1polity2*Wmember+ l1log.pop + l1log.gdppc  +l1bpc+l1ecris+l1yrsoffc +l1gatt.wto
                  + l1usheg + l1avetar + Eight + wtrade+(1|cowcode) + (l1polity2|year), data=traderep)

trade.mlm <- lmer(newtar~l1polity2*Wmember+ l1log.pop + l1log.gdppc  +l1bpc+l1ecris+l1yrsoffc+l1gatt.wto
                  + l1usheg + l1avetar + Eight + wtrade+(1+Wmember|l1polity2)+ (1|cowcode) + (l1polity2|year), data=traderep)

trade.mlm2 <- lmer(newtar~l1polity2*Wmember+ l1polity2*l1usheg + l1polity2*l1avetar + l1polity2*Eight +
                   l1polity2*wtrade++ l1log.pop + l1log.gdppc  +l1bpc+l1ecris+l1yrsoffc+l1gatt.wto+ 
                   (1|cowcode) + (l1polity2|year), data=traderep)

trade.mlm3 <- lmer(newtar~l1polity2*Wmember + l1polity2*l1avetar+ l1log.pop + l1log.gdppc  +l1bpc+l1ecris+l1yrsoffc+l1gatt.wto
                  + l1usheg + Eight + wtrade+ (1|cowcode) + (l1polity2|year), data=traderep)

##################################
# run Type IV with heterogeneity with the data 1970-2008
##################################
  load("trade011312.RData")
  traderep <- trade011312
  load("WmWTOlist.RData")
  # the previous one
  # traderep <- read.csv("tradeupdateimp2.csv")
  # W <- read.csv("distanceW.csv")
  attach(traderep)
  source("functions.r")
  source("BSMLMWT.R")
  Y <- newtar
  n <- length(Y)
  X1 <- cbind(1, l1polity2, l1log.pop, l1log.gdppc,l1bpc,l1ecris,l1yrsoffc,
              l1gatt.wto,l1usheg, l1avetar, Eight, Wmember, wtrade)
  Si <- rep(1, n)
  Dt <- rep(1, n)
  Uid <- cowcode
  Tid <- traderep$year
  W <- WmWTOlist
  WT <- 1:39
  rho.start <- rep(0, 39)
  precision <- 1
  beta0 <- rep(0, 13)
  B0 <- diag(400, 13)
  a0 <- 1
  b0 <- 1
  E0 <- 1
  e0 <- 1
  D0 <- 1
  d0 <- 1
  mcmc <- 10000
  burnin <- 5000
  thin <- 1
  tracking <- 50
  model3 <-   BSMLMWT.MCMC(Uid,Tid, Y, X1, Si, Dt, Ai="NULL", Zt="NULL",W, WT, timeint.add=FALSE, unitint.add=FALSE,
                       rho=rho.start, precision, beta0, B0, a0, b0,  D0, d0, E0, e0,
                       mcmc, burnin, thin, tracking)
  NewSpatialMCMC <- model3
  save(NewSpatialMCMC, file="NewSpatialMCMC.RData")


#########################
# Run Type V interaction with data from 1970-2008
#########################
 load("trade011312.RData")
  traderep <- trade011312
  load("WmWTOlist.RData")
  # the previous one
  #traderep <- read.csv("tradeupdateimp2.csv")
  # W <- read.csv("distanceW.csv")
  attach(traderep)
  source("functions.r")
  source("BSMLMWT.R")
# [1] "l1log.pop"   "l1log.gdppc" "l1bpc"       "l1ecris"     "l1yrsoffc"  
# [6] "l1gatt.wto"  "1"       l1usheg, l1avetar, Eight, Wmember, wtrade
# l1polity2, interactions
  Y <- newtar
  n <- length(Y)
  X1 <- cbind(l1log.pop, l1log.gdppc,l1bpc,l1ecris,l1yrsoffc,
              l1gatt.wto)
  Si <- rep(1, n)
  Dt <- cbind(1, l1polity2)
  Zt <- cbind(1, l1usheg, l1avetar, Eight, Wmember, wtrade)
  Uid <- cowcode
  Tid <- traderep$year
  W <- WmWTOlist
  WT <- 1:39
  rho.start <- rep(0, 39)
  precision <- 1

  beta0 <- rep(0, 18)
  B0 <- diag(400, 18)
  a0 <- 1
  b0 <- 1
  E0 <- diag(2, 2)
  e0 <- 3
  D0 <- 1
  d0 <- 1
  mcmc <- 10000
  burnin <- 5000
  thin <- 1
  tracking <- 50

  model2 <-   BSMLMWT.MCMC(Uid,Tid, Y, X1, Si, Dt, Ai="NULL", Zt,W, WT, timeint.add=FALSE, unitint.add=FALSE,
                       rho=rho.start, precision, beta0, B0, a0, b0,  D0, d0, E0, e0,
                       mcmc, burnin, thin, tracking)
  SpatialMLMMCMC <- model2
  save(SpatialMLMMCMC, file="SpatialMLMMCMC.RData")

#\\\\\\\\\\\\\\\\\\\\\\\ Code Ends Here \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\













